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This  paper  provides  the  first  analysis  of  the  seasonal  evolution  of  the  sonic  layer  depth  (SLD)  relative  to  the 
mixed  layer  depth  (MLD)  for  the  Aegean,  Marmara,  Black,  and  Azov  Seas.  SLD  identifies  the  acoustic  ducting 
capabilities  of  the  upper  ocean  and  is  of  interest  to  investigations  of  upper  ocean  acoustics.  A  monthly  SLD 
climatology  on  a  regular  0.25°  x  0.25°  grid  is  constructed  from  interpolation  of  available  quality-controlled 
ocean  temperature  and  salinity  profiles  using  the  kriging  methodology.  A  four  step  pre-processing  procedure 
is  designed  to  reduce  noise  and  the  effects  of  sampling  irregularities.  Monthly  SLD  fields  are  then  compared 
relative  to  the  much  more  widely  studied  MLD  as  computed  using  four  different  methods  from  the  recent 
scientific  literature.  The  goals  of  this  analysis  are  to  characterize  the  SLD  relative  to  the  MLD  and  provide  a 
means  for  computing  SLD  from  limited  hydrographic  information  and/or  MLD  estimates.  Very  deep  SLD 
values  are  found  during  winter,  in  the  Aegean  and  Black  Seas,  when  the  near  surface  temperature  values 
become  lower  than  the  temperature  below  the  permanent  pycnocline.  When  this  occurs,  the  SLD  drops  to 
the  bottom  while  the  MLD  remains  much  shallower  at  the  seasonal  pycnocline.  For  the  months  of  May 
through  October  the  SLD  tends  to  be  less  than  25  m  for  the  entire  region.  It  is  demonstrated  that  MLD 
obtained  from  the  four  methodologies  have  high  correlations  with  SLD  over  the  annual  cycle,  indicating  a 
robust  relationship.  As  a  result,  SLD  can  be  estimated  using  least  squares  regression  coefficients  when  salinity 
is  unavailable  or  when  observation  profiles  do  not  extend  to  deeper  levels. 
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1.  Introduction 

Conditions  of  the  upper  ocean  often  result  in  a  near  surface 
acoustic  duct  that  limits  the  downward  transmission  of  sound  and 
results  in  the  acoustic  spreading  to  be  approximately  cylindrical 
(Urick,  1983).  Sound  trapped  in  a  surface  duct  is  primarily  transmitted 
outward  from  the  source  in  an  expanding  disk.  This  occurs  when  the 
near  surface  environment  is  upward  refracting  due  to  an  increase  of 
sound  speed  with  depth.  The  depth  over  which  the  increasing  sound 
speed  penetrates  is  called  the  sonic  layer  depth  (SLD)  (e.g.,  Kerman, 
1993).  The  acoustic  frequency  that  can  propagate  in  this  type  of  duct  is 
dependent  on  the  SLD.  For  a  given  SLD  there  is  a  cutoff  frequency 
above  which  sound  will  be  trapped  in  the  surface  duct.  Sound  with  a 
lower  frequency  is  not  trapped  and  will  spread  spherically  in  all 
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directions,  leading  to  more  rapid  horizontal  attenuation  (Etter,  2003). 
Sound  in  a  duct  (with  cylindrical  spreading)  can  be  transmitted  much 
farther  horizontally  than  in  non-ducting  environments  (with  sphe¬ 
rical  spreading).  This  phenomenon  substantially  influences,  for 
example,  ocean  communications  (Siderius  et  al.,  2007),  acoustic 
tomography  (Sutton  et  al.,  1993),  and  naval  operations  (Urick,  1983). 

A  parameter  related  to  SLD,  the  mixed  layer  depth  (MLD),  is  also 
important  because  it  identifies  the  penetration  depth  of  turbulence 
near  the  ocean  surface  and  has  a  wide  array  of  influence  on  a  variety  of 
upper  ocean  processes  from  air-sea  exchange  (e.g.,  Chen  et  al.,  1994) 
to  biological  interactions  (e.g.,  Siegel  et  al.,  2002).  For  a  typical  water 
column,  where  isothermal  and  isohaline  surface  layer  depths  are  equal 
and  temperature  decreases  below  that  depth,  the  MLD  and  SLD  are 
equal.  This  is  because  sound  speed  increases  with  depth  due  only  to 
increasing  pressure  until  a  local  maximum  at  the  MLD,  below  which 
sound  speed  decreases  with  decreasing  temperature.  There  are  other 
conditions,  however,  where  the  SLD  may  be  lower  or  higher  in  the 
water  column  than  the  MLD.  The  characteristics  of  how  SLD  and  MLD 
become  different  are  discussed  in  this  paper  because  they  provide 
useful  information  regarding  the  structure  of  the  upper  ocean 
temperature  and  salinity  fields. 

In  a  recent  analysis  (Helber  et  al.,  2008),  the  global  deviations  of  the 
MLD  relative  to  the  SLD  were  investigated  for  the  annual  cycle.  It  was 
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found  that  in  the  spring  when  new  stratification  events  occur,  the  SLD 
can  be  substantially  deeper  than  the  MLD.  In  the  present  study,  we 
examine  differences  between  SLD  and  MLD  in  more  detail  in  a  region 
including  a  small  part  of  the  eastern  Mediterranean  Sea  and  the 
Aegean,  Marmara,  Black  and  Azov  Seas.  This  is  accomplished  by 
interpolating  the  SLD  and  MLD  values  to  a  0.25  x  0.25  grid  using  kriging 
(Diggle  and  Ribeiro,  2007).  This  region  is  chosen  because  there  have 
not  been  any  previous  studies  quantifying  and  examining  SLD  features. 

The  Black  Sea  is  heavily  influenced  by  river  input  that  along  with 
precipitation  exceeds  evaporation  losses  and  is  balanced  by  outflow 
across  the  Bosphorus  (Ozsoy  and  Unliiata,  1998).  As  a  result,  the  Black 
Sea  is  characterized  by  fresher  surface  water  residing  above  warmer 
saltier  deep  water  (a  condition  associated  with  “diffusive  convection,” 
e.g.,  Kantha  and  Clayson,  2000).  Above  a  permanent  pycnocline  there 
exists  a  cold  intermediate  layer  (CIL)  that  is  replenished  in  the  winter 
when  the  seasonal  thermocline  is  deepest  (Ozsoy  and  Onliiata,  1998). 
The  Aegean  Sea  is  influenced  by  cooler  fresher  surface  water 
originating  in  the  Black  Sea  that  passes  through  the  Bosphorus  and 
the  Marmara  Sea  and  enters  the  Aegean  through  the  Dardanelles 
Strait.  This  water  flows  across  the  northern  Aegean  shelf  and 
southward  along  the  western  coastline  (Kourafalou  and  Barbopoulos, 
2003;  Olson  et  al.,  2007).  Cold,  saltier  Mediterranean  water  enters  the 
lower  layers  of  the  Marmara  Sea  through  the  Dardanelles  Strait 
(Be^iktepe  et  al.,  1994).  The  Sea  of  Azov  is  very  shallow,  with  an 
average  depth  of  5  m,  but  stratification  does  exist  (Debolskaya  et  al., 
2008).  In  all  of  these  regions,  to  some  degree,  the  surface  layers 
defined  by  MLD  and  SLD  will  be  associated  with  cooler,  fresher  surface 
water  above  saltier,  warmer  subsurface  water.  As  will  be  shown,  the 
evaluation  of  the  SLD  relative  to  the  MLD  is  impacted  by  this  unique 
hydrographic  structure,  making  this  an  excellent  region  for  studying 
SLD  and  its  relation  to  MLD. 

Estimating  SLD  from  observation  profiles  is  achieved  using  a 
method  derived  by  Helber  et  al.  (2008).  We  use  four  methodologies 
for  MLD  taken  from  the  recent  literature  (Kara  et  al.,  2000b;  Lorbacher 
et  al.,  2006)  that  are  discussed  by  Kara  et  al.  (2009).  For  two  of  the 
algorithms,  temperature-only  profiles  are  used,  while  the  other  two 
select  the  MLD  from  density  profiles.  Comparison  of  SLD  with  both 
temperature  and  density  MLD  estimates,  all  interpolated  on  a  0.25° 
grid,  provide  information  regarding  the  salinity  contribution  to  the 
upper  ocean  in  addition  to  characterizing  upper  ocean  acoustics  and 
the  penetration  depth  of  surface  turbulence.  Use  of  several  MLD  meth¬ 
ods  enhances  the  results  by  ensuring  that  the  observed  phenomena  are 
not  an  artifact  of  one  MLD  methodology. 

The  differences  between  the  SLD  and  MLD  arise  because  of 
differences  in  the  sensitivities  of  density  and  sound  speed  to 
temperature,  salinity,  and  pressure.  For  density,  salinity  variability  is 
often  small  enough  that  the  MLD  is  controlled  by  temperature. 
Salinity,  however,  can  have  a  large  impact  producing,  for  example,  the 
“barrier  layer”  (e.g.  Kara  et  al.,  2000a),  where  the  salinity  contribution 
results  in  a  MLD  that  is  shallower  than  temperature  alone  would 
dictate.  Salinity  has  a  smaller  impact  on  sound  speed  relative  to 
density  because  sound  speed  is  more  sensitive  to  temperature  than 
salinity.  A  fundamental  difference  between  sound  speed  and  density, 
and  therefore  SLD  and  MLD,  is  that  sound  speed  does  not  influence 
turbulence  whereas  density  does.  These  differences  results  in  a 
relationship  between  SLD  and  MLD  that  we  explore  in  this  analysis. 

A  major  goal  of  this  paper  is  to  identify  the  conditions  under  which 
SLD  and  MLD  differ  in  time  and  space  in  the  Aegean,  Marmara,  Black 
and  Azov  Seas  region.  In  addition,  statistical  relationship  between  SLD 
and  MLD  is  exploited  in  order  to  provide  a  means  for  predicting  SLD 
from  MLD  estimates.  This  capability  is  particularly  useful  when  only 
temperature  profiles  or  MLD  forecasts  are  available. 

In  Section  2  we  describe  the  sources  of  the  observation  profiles  and 
the  quality  control  and  data  editing  procedures.  The  methods  for 
estimating  MLD,  SLD  and  computing  the  interpolated  fields  are 
described  in  Section  3.  The  monthly  SLD  climatology  is  described  in 


Section  4  and  the  differences  between  SLD  and  MLD  are  discussed  in 
Section  5.  Section  6  describes  the  statistical  relationship  of  SLD  and 
MLD  and  Section  7  concludes  the  paper.  The  Appendix  A  contains 
information  about  the  data  preparation  steps  taken  for  producing  the 
interpolated  SLD/MLD  0.25°  gridded  fields. 

2.  Data 

The  available  historically  observed  temperature  profiles  (T-only) 
and  paired  profiles  of  temperature  and  salinity  (T  and  S)  are  acquired 
from  the  following  sources:  (1)  the  World  Ocean  Database  2005 
(WOD05)  (Boyer  et  al.,  2006),  (2)  the  U.  S.  Navy's  Master  Oceano¬ 
graphic  Observation  Data  Set  (MOODS)  (Teague  et  al.,  1990),  (3)  Argo 
(Gould  et  al.,  2004),  and  (4)  the  National  Oceanographic  Data  Center 
(for  the  Sea  of  Azov,  Matishov  et  al.,  2006).  The  reader  is  referred  to 
Table  1  for  the  number  of  profiles  of  T-only  and  T  and  S  from  each 
source.  The  first  three  sources  of  data  are  partially  redundant  but  all 
are  considered,  because  after  the  removal  of  duplicates  there  are 
unique  data  in  each.  The  WOD05  has  less  stringent  quality  control 
(QC)  procedures  than  MOODS  but  potentially  more  data  sources. 
There  is  an  approximately  6%  increase  in  data  volume  with  the 
addition  of  WOD05  data  relative  to  the  MOODS  data  set  alone.  Older 
Argo  data  are  included  in  the  MOODS  and  WOD05,  and  the  latest  Argo 
data  from  the  Data  Acquisition  Centers  (DAC)  have  been  included  with 
the  application  of  the  latest  QC  recommendations.  As  an  example  of 
the  data  coverage  for  observations  with  profile  pairs  of  T  and  S,  Fig.  1 
shows  the  locations  of  usable  profiles  for  all  years  in  the  month  of 
October  only.  Note  that  for  the  Marmara  and  Azov  Seas  there  are  very 
few  data  points,  and  therefore  the  results  within  these  seas  is  less 
robust  than  in  the  Black  and  Aegean  Seas. 

While  the  data  are  quality-controlled,  there  are  still  some  existing 
errors  that  include  profile  location,  XBT  drop  rate,  low  vertical 
resolution,  etc.  In  the  interest  of  including  as  much  data  as  possible, 
profiles  with  relatively  low  vertical  resolution  or  large  gaps  in  the 
vertical  have  been  retained  in  the  data  set.  To  help  minimize  these  and 
other  sources  of  potential  errors,  SLD  “outliers”  are  removed  statis¬ 
tically  as  explained  in  Appendix  A.  Errors  that  are  not  also  “outliers” 
are  statistically  undetectable  and  therefore  cannot  be  removed  in  this 
manner. 

3.  Methods 

The  SLD  algorithm  is  described  in  detail  in  Helber  et  al.  (2008)  and 
will  be  explained  briefly  here.  The  first  step  in  estimating  the  SLD  from 
T  and  S  profile  pairs  is  the  computation  of  sound  speed  at  each  depth 
level  using  the  nine-term  equation  of  Mackenzie  (1981): 

c  =  1448.96  +  4.591T  -  5.304  x  +  2.374  x 

+  1.340(S  -  35)  +  1.630  x  10“^D  +  1.675  x  10“’D^ 

-  1.025  x  10“^T(S  -  35)  -  7.139  x  10“’^TD^ 

Sound  speed,  c  (m  s“’),  is  a  function  of  temperature,  T  (°C), 
salinity,  S  (psu),  and  depth,  D  (m)  and  is  non-linear.  The  first  three 
non-constant  terms  depend  only  on  T  and  are  the  largest  in  the  upper 
ocean.  Only  at  depth  do  the  terms  including  D  have  large 

Table  1 

The  number  of  total,  T-oniy,  and  T  and  S  profiles  from  each  source  in  the  analysis  domain 
of  Fig.  1  (22.00  E  to  41.84  E  and  34.50310  N  to  47.31377  N). 
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Fig.  1.  Observation  locations  for  T  and  5  profiles  with  at  least  three  usable  depth  levels  over  the  study  region  in  the  month  of  October.  Profiles  with  7-only  are  more  numerous  and  are 
not  shown. 


contributions.  As  noted  earlier,  sound  speed  is  less  sensitive  to  5  than 
to  r,  as  compared  to  density  for  typical  oceanic  values  for  T  and  5. 

Starting  at  the  surface  downward,  each  depth  level  {k  =  1,2,3,-  K)  of 
a  given  profile  is  examined  using  two  criteria: 

( 1 )  Is  the  kth  sound  speed  value  greater  than  those  at  levels  less 
than  k 

(2)  Is  the  kth  sound  speed  value  greater  than  that  at  level  k  + 1 

Each  depth  level  that  satisfies  the  above  2  criteria  identities  a  local 
depth  maximum  in  sound  speed  and  an  overall  upward  refracting 
(negative)  sound  speed  gradient  from  the  surface  to  the  kth  depth  level. 
If  more  than  one  data  level  meets  these  criteria,  several  other  tests  are 
conducted  to  identify  the  nature  of  the  sound  speed  gradient  above  each 
level.  In  most  cases  the  deepest  identified  level  that  meets  the  criteria  is 
taken  as  the  SLD,  but  there  are  additional  tests  conducted  to  eliminate 
potential  SLD  errors  that  arise  with  these  simple  criteria  and  to  omit 
solutions  that  are  inappropriate  for  a  particular  application. 

The  tests  are  employed  using  two  tuning  parameters  that  rep¬ 
resent  the  average  frequency  (on  a  log  scale)  and  the  lowest  frequency 
expected  for  an  acoustic  application.  For  this  study  we  chose  a  low 
frequency  range  generally  used  in  Navy  operations  (Kunz  and  Gauss, 
2000)  defined  by  the  mid  and  lowest  frequencies  of  680  Hz  and  70  Hz, 
respectively.  For  example,  a  test  is  conducted  to  check  for  the  exis¬ 
tence  of  a  subsurface  duct  above  each  of  the  identified  data  levels. 
In  cases  where  there  is  a  subsurface  duct,  the  SLD  is  non-existent  and 
is  set  to  zero.  To  determine  the  presence  of  a  subsurface  duct,  the 
minimum  cutoff  frequency  is  computed  by  integrating  the  sound 
speed  profile  (Etter,  2003).  If  the  frequency  is  greater  than  the  mid 
frequency  (680  Hz)  then  the  duct  above  the  identified  level  is  large 
enough  to  negate  the  existence  of  a  surface  duct  and  SLD  is  set  to  zero. 
Additional  tests,  to  ensure  the  proper  SLD  is  identified,  are  described 
by  Helber  et  al.  (2008). 

Four  different  MLD  methodologies  are  employed  to  help  identify 
only  the  most  robust  relationships  between  SLD  and  MLD.  We  have 
chosen  two  algorithms  that  are  based  solely  on  T-only  profiles  and 
two  algorithms  that  require  T and  S  profiles.  As  described  by  Kara  et  al. 
(2009)  the  T-only  methods  are  an  isothermal  threshold  method  (Kara 
et  al.,  2000b)  denoted  MLD  (I)  and  a  curvature  based  method 


(Lorbacher  et  al.,  2006)  denoted  MLD(L).  The  methods  that  require 
both  T and  S  to  compute  at  are  a  variable  Ot  threshold  method  (Kara  et 
al.,  2000b)  denoted  MLD(ctt)  and  a  constant  density  threshold 
method  denoted  MLD(cJc)  employed  using  the  programs  of  Kara  et 
al.  (2000b).  Hereafter,  MLD  will  be  used  to  represent  the  mixed  layer 
depth  in  general,  while  the  notation  MLD(L),  MLD(r),  MLD(crT),  and 
MLD(o'c)  will  represent  the  MLD  as  estimated  from  observations  using 
the  associated  methodologies  listed  above. 

To  create  monthly  gridded  fields  of  SLD  and  MLD  we  have  processed 
the  data  to  reduce  noise  and  normalized  the  error  levels  using  several 
steps  described  in  the  Appendix  A.  These  climatologically  processed 
observations  represent  an  unknown  underlying  process  and  contain 
unavoidable  noise.  The  noise  comes  from,  for  example,  interannual 
variability,  systematic  and  random  errors,  and  inhomogeneous  sam¬ 
pling.  With  this  in  mind,  interpolation  of  the  irregular  data  locations  are 
then  performed  using  kriging.  Kriging  takes  into  account  the  presence  of 
noise  and  constructs  a  linear  predictor  based  on  the  covariance  structure 
of  the  data  (Cressie,  1991;  Diggle  and  Ribeiro,  2007).  The  kriging 
parameters  that  control  the  linear  predictor  and  the  amount  of 
smoothing  applied  to  the  output  field  are  the  range  cj),  sill  and 
nugget.  Using  an  empirical  technique,  the  sill  and  nugget  are  linearly  fit 
to  the  biweight  standard  deviation  of  the  data  using  a  range  of  2.5°.  The 
nugget  is  chosen  such  that  the  amount  of  smoothing  applied  has  a 
characteristic  length  equal  to  1.4°.  The  details  of  this  procedure  are 
described  in  Kara  et  al.  (2009). 

4.  Monthly  SLD  climatology 

Following  methods  in  the  Appendix  A  and  the  kriging  interpolation 
described  in  the  Section  3  and  in  Kara  et  al.  (2009),  we  now  present  the 
monthly  mean  SLD  fields  and  compare  them  to  the  four  MLD  fields 
mentioned  above.  It  can  be  seen  in  Fig.  2  that  the  SLD  has  a  seasonal 
cycle  with  relatively  deep  SLD  during  months  of  January,  February,  and 
March.  February  has  the  greatest  extent  of  very  deep  SLDs.  In  the 
southern  part  of  the  Aegean  Sea  and  the  south  western  part  of  the  Black 
Sea  the  SLD  is  deeper  than  200  m  during  January  and  February.  During 
March  the  deep  SLDs  are  restricted  to  the  center  of  the  western  Black 
Sea  gyre  and  isolated  regions  of  the  Aegean  Sea.  Because  of  the  deep 
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SLD  values,  January,  February,  and  March  would  be  good  months  for 
long  range  upper  ocean  acoustic  communication. 

The  very  deep  SLD  values  occur  at  these  locations  during  winter 
when  the  surface  temperatures  are  cooler  than  the  deep  water 
temperatures.  There  is  a  permanent  pycnocline  at  50-100  m  in  the 
Black  Sea  that  occurs  at  the  bottom  of  the  CIL  (Ozsoy  and  Unliiata, 
1998)  and  below  500  m  the  water  of  Mediterranean  origin  is 
essentially  well-mixed  (Poulain  et  al.,  2005).  When  the  near  surface 
temperatures  are  below  approximately  8  °C  (or  the  temperature  of  the 
CIL),  the  sound  speed  increases  with  depth  all  the  way  to  the  bottom 
and  consequently  the  SLD  is  also  at  the  bottom  (Fig.  2).  An  example 
profile  is  shown  in  Fig.  3. 

This  phenomenon  has  sharp  transitions  in  space  and  time.  As  near 
surface  temperatures  seasonally  become  cooler  in  the  fall,  the  SLD 
becomes  gradually  deeper  until  the  surface  temperatures  dip  below 
the  temperatures  of  the  CIL.  At  this  point  the  SLD  immediately  drops 
to  the  bottom.  This  shift  to  deep  SLDs  occurs  in  short  distances 
horizontally  as  can  be  seen  in  Fig.  2  for  months  of  January,  February, 
and  March  where  large  SLD  gradients  occur  around  the  deep  SLD  (red) 
areas.  The  standard  deviation  of  SLD  is  also  large  during  the  months 
where  deep  SLDs  occur  (Table  2). 


For  these  deep  SLD  cases  there  is  some  ambiguity.  Many  CTD  casts 
do  not  go  to  the  bottom.  For  a  sound  speed  profile  that  is  upward 
refracting  over  the  whole  water  column,  the  actual  SLD  would  be  at 
the  ocean  bottom  but  the  analysis  of  the  observation  would  put  the 
SLD  at  the  bottom  of  the  CTD  cast.  For  this  reason.  Fig.  2  does  not  show 
the  SLD  variability  when  the  SLD  is  deeper  than  200  m. 

The  seasonal  cycle  of  SLD  can  also  be  seen  in  the  time  series  at 
various  point  locations  (Fig.  4).  The  deepest  SLDs  and  the  largest 
seasonal  variability  occur  in  the  western  gyre  of  the  Black  Sea  in 
February  (see  also  Table  2).  The  smallest  variability  is  during  the 
months  of  May  through  September  and  in  the  relatively  shallow  Sea  of 
Azov  (Fig.  2).  The  Aegean  Sea  has  a  longer  deep  SLD  season  that  starts 
sooner  and  ends  later  than  in  the  Black  Sea. 

5.  Differences  between  SLD  and  MLD 

Comparisons  between  different  types  of  surface  layers  provide 
information  regarding  the  upper  ocean  structure.  For  example,  when 
the  SLD  is  deeper  than  the  MLD,  salinity  and  temperature  both  tend  to 
increase  with  depth  (Fig.  3).  Salinity  must  increase  with  depth  at  a 
sufficient  rate,  because  otherwise  the  profile  would  be  unstable.  When 


Fig.  2.  Monthly  mean  SLD  climatology  in  the  Aegean,  Marmara,  Black  and  Azov  Seas  by  month.  SLD  fields  are  produced  based  on  individual  T  and  S  profiles  available  from  WOD05, 
MOODS,  and  Argo  data  archives.  See  text  for  further  details.  The  color  palette  is  the  same  for  all  months  to  show  seasonal  variability. 
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Fig.  3.  An  example  profile  where  the  SLD  exists  at  the  bottom  of  the  profile,  SLD  =  522  m. 
The  location  of  the  profile  is  30°N,  43°E  taken  on  19/02/1975  and  MLD(r)  =  141,  MLD(aT)  = 
53,  MLD(ac)  =  55,  and  MLD(L)  =  10  m. 


the  MLD  is  deeper  than  the  SLD,  it  is  because  SLD  responds  to  small 
sound  speed  maxima  in  profiles  that  are  relatively  uniform  with  depth 
(Fig.  6).  To  explore  these  relationships  further,  we  investigate  whether 
or  not  SLD  is  typically  deeper  than,  shallower  than,  or  identical  to  MLD. 

Monthly  mean  MLD  climatologies  based  on  the  four  MLD 
methodologies  discussed  in  Section  3  are  used  to  compute  differences 
between  SLD  and  MLD.  As  evident  from  Fig.  5,  SLD  can  be  deeper  or 
shallower  than  MLD  depending  on  the  month  of  the  year.  This  is  true 
regardless  of  which  method  is  used  for  computing  MLD.  There  are  two 
main  features  in  the  SLD  versus  MLD  differences.  The  SLD  is  typically 
much  deeper  than  MLD  in  the  southwestern  Black  Sea  and  parts  of 
the  Aegean  Sea  in  January  (Fig.  5a)  and  February  and  March  (not 
shown).  When  SLD  drops  to  the  bottom,  as  described  in  Section  4, 
the  MLD  remains  with  the  seasonal  pycnocline.  During  February 
the  difference  is  so  large  and  widespread  that  it  even  dominates  the 
basin-wide  average  (Table  3).  This  difference  is  enhanced  for  the 
MLD(L)  field  due  to  its  wintertime  shallow  bias  (Kara  et  al.,  2009). 
The  difference  in  SLD-MLD(L)  is  on  average  77  m  over  the  study 
region  (Table  3).  The  variance  of  the  basin-wide  difference  fields  is 
considerably  larger  during  winter  and  early  spring  due  to  the  vari¬ 
ance  in  SLD  (Table  3). 

The  other  main  feature  is  that  MLD  is  deeper  than  SLD  over  the 
entire  Black  Sea  for  the  rest  of  the  year,  but  the  magnitude  of  this  bias 
is  generally  smaller.  The  values  for  the  difference  SLD-MLD  in  Table  3 
are  mostly  negative  for  all  months  except  February.  The  MLD(L) 
results  are  the  exception.  An  area  where  MLD  can  be  considerably 
deeper  than  SLD  is  in  the  spring  when  new  stratification  occurs  in  the 
Aegean  Sea  (Fig.  5b).  This  bias  is  absent  for  MLD(L)  because  the 
curvature  based  methodology  detects  small  temperature  deviations 
associated  with  new  stratification  (Helber  et  al.,  2008).  These  small 


Table  2 

Monthly  standard  deviation  of  the  SLD  for  the  Black  Sea  between  latitudes  41  °N  and 
45°N. 


Month 

Jan 

Eeb 

Mar 

Apr 

May 

Jun 

Jui 

Aug 

Sep 

Oct 

Nov 

Dec 

SLD  standard 
deviation,  m 

38 

212 

128 

23 

7 

4 

6 

9 

13 

16 

23 

29 

Outlier  removal  using  the  Z-score  statistic  (see  Appendix  A)  was  applied  prior  to  the 
computation  of  the  standard  deviation. 


Fig.  4.  Climatological  mean  time  series  of  SLD  at  three  locations  over  the  annual  cycle. 
These  locations  are  randomly  chosen  to  represent  typical  SLD  values  in  the  Aegean 
(39°N  24°E),  Black  (43°N  30°E)  and  Azov  Seas  (46°N  37°E). 

deviations  are  also  sufficient  to  create  a  SLD  at  that  depth  (Fig.  6). 
Since  the  variance  of  temperature  with  depth  is  relatively  small,  the 
MLD((Jt)  and  MLD((Jc)  values  remain  deeper.  The  MLD(L)  values 
generally  remain  closer  to  SLD  than  the  MLD((Jt)  and  MLD((Jc)  values 
except  for  during  the  deep  winter  SLD  regions.  During  the  boreal 
summer  and  fall,  SLD  remains  shallower  than  all  four  MLDs  (Fig.  5c,d). 

To  illustrate  the  differences  between  SLD  and  MLD((Jt)  in  another 
way,  we  produce  a  scatter  diagram  in  January  (Fig.  7)  and  July  (Fig.  8). 
The  kriging  interpolation  field  data  are  used  to  construct  the  figures. 
For  January  there  exists  a  statistically  significant  linear  correlation 
coefficient  of  0.74.  There  are  two  off  diagonal  “lobes:”  (1)  in  the  left 
side  of  Fig.  7  where  MLD(o't)  is  between  50  and  100  m  that  represents 
values  primarily  in  the  Black  Sea  (Fig.  5)  and  (2)  the  lower  right  side  of 
Fig.  7  where  SLD  is  between  50  and  125  m  that  represents  values  in  the 
Mediterranean. 

For  July  there  is  no  statistically  significant  linear  correlation, 
because  there  are  two  competing  regions  with  different  character¬ 
istics.  One  region  is  where  SLD  is  essentially  zero  and  the  other  is 
where  SLD  is  non-zero  but  smaller  than  MLD((Jt).  When  SLD  is  not 
zero,  it  is  still  shallower  than  MLD(o't).  Regionally,  correlation 
between  SLD  and  MLD  is  significantly  higher,  as  will  be  shown  in 
Section  6.  This  is  evident  in  the  scatter  plots  of  Figs.  7  and  8  because 
each  of  the  “lobes”  represents  a  region  of  high  correlation. 

6.  Relationship  between  SLD  and  MLD 

In  this  section  we  seek  a  simple  relationship  between  SLD  and  MLD 
values  over  the  seasonal  cycle.  This  could  help  identify  regions  where 
SLD  can  be  obtained  from  T-only  profiles  using  MLD(L)  or  MLD  (I)  and 
a  simple  computation,  e.g.,  adding  or  subtracting  a  constant  offset. 
This  could  be  advantageous  over  the  alternative  of  using  T  and  S 
profiles,  which  are  generally  less  available  than  T-only  profiles.  In 
addition,  for  T-only  observations  that  extend  less  than  200  m,  as  for 
many  bathythermographs,  such  a  relationship  could  be  used  to 
estimate  the  SLD  that  could  be  much  deeper.  We  use  all  MLD 
climatologies  in  seeking  a  relationship  to  test  the  sensitivity  of  the 
statistics  to  the  MLD  definition. 

We  evaluate  time  series  of  SLD  and  MLD  over  the  seasonal  cycle. 
The  statistical  relationships  used  in  comparisons  between  the 
12  monthly  MLD  (X)  and  SLD  (Y)  values  at  each  grid  point  are 
expressed  as  follows: 

(^i  - / (w),  (1) 

"i=l 

where 

Y  =  a  +  bX  +  c  (2) 

and  n  =  12,  R  is  the  correlation  coefficient,  and  X  (Y)  and  dx  (oy)  are 
the  mean  and  standard  deviations  of  MLD  (SLD)  values,  respectively. 


S186 


R.W.  Helber  et  al.  /  Journal  of  Marine  Systems  78  (2009)  S181-S190 


(a)  Jan 


(b)Apr  (c)Jul 


(d)  Oct 


Fig.  5.  Differences  between  SLD  and  the  four  MLD  definitions  MLD(aT),  MLD(ac),  IVILD(T),  and  IVILD(L).  Results  are  shown  in  (a)  January,  (b)  April,  (c)  July  and  (d)  October,  separately. 
Note  the  color  palette  for  the  difference  fields  in  January  and  April  is  different  from  that  in  July  and  October. 


In  the  regression  Eq.  (2),  Y  is  the  dependent  variable,  X  is  the  inde¬ 
pendent  variable  (or  covariate),  a  is  the  intercept,  b  is  the  slope  or 
regression  coefficient,  and  c  is  the  error  term.  The  regression  equation 
will  specify  the  average  magnitude  of  the  expected  change  in  SLD  for 
the  given  MLD  climatology. 


We  first  determine  how  SLD  and  MLD  are  correlated  over  the 
seasonal  cycle.  The  strength  of  the  linear  relationship  between  SLD 
and  MLD  is  determined  by  the  R  value  (Eq.  (1)),  which  ranges  from 
—  1  to  1.  An  value  of  —1  (1)  indicates  a  very  strong,  negative 
(positive)  linear  relationship,  and  an  R  value  of  0  indicates  no  linear 


Table  3 

SLD  differences. 


Basin  averaged  differences  of  various  mixed  layer  depths  subtracted  from  SLD. 

Standard  deviations  are  calculated  from  the  entire  domain  (22.00  to  41.84  E,  34.5031  to  47.31377  N)  of  kriged  sea  points. 
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Fig.  6.  An  example  profile  during  April  1990  where  the  MLD(L)  (53  m)  and  the  SLD  (53  m)  are  relatively  shallow  while  MLD(aT)  (201  m),  MLD(ac)  (201  m)  and  MLD(T)  (201  m)  are 
considerably  deeper. 


relationship.  The  resulting  correlations  given  in  Fig.  9  clearly  reveal 
strong  and  positive  values  (i^>0.6)  over  most  of  the  domain.  This  is 
true  regardless  which  MLD  methodology  is  used. 

To  decide  whether  or  not  a  correlation  value  represents  a  good  or 
poor  linear  relationship,  we  apply  the  student's  t-test.  Student's  t-test 
is  applied  to  test  the  null  hypothesis  that  two  variables  (i.e.,  SLD  and 
MLD)  are  not  correlated  at  a  given  grid  point.  For  the  n  =  12  sample 
values  in  this  paper,  the  bivariate  normal  distribution  of 
tn-2  =  RVn  —  2  /  must  have  an  absolute  value  of  at  least 

0.53  for  R  to  be  statistically  significant  at  the  95%  {a  =  0.05 )  confidence 
level. 

Given  the  significant  correlation  values  over  almost  the  entire 
domain,  including  the  Aegean,  Marmara,  Black  and  Azov  Sea,  a  linear 
regression  can  be  used  to  examine  the  relationship  between  SLD  and 
MLD.  Details  of  this  commonly-used  linear  regression  technique  may 
be  found  in  Wilks  (1995).  Our  objective  is  to  fit  a  straight  line  to  values 
between  SLD  and  MLD  obtained  at  each  0.25°  x  0.25°  over  the  seasonal 
cycle.  We  would  like  to  have  a  line  that  is  in  some  sense  closest  to  all  of 


the  data  points  simultaneously,  i.e.,  a  method  of  finding  a  line  via  the 
least  squares  approach.  This  process  is  done  at  each  grid  point  based 
on  12-month  SLD  and  MLD  time  series.  Slope  and  intercept  values  of 
the  least  squares  line  are  then  computed  as  described  in  Eq.  (2). 

Fig.  9  shows  that  slope  values  are  remarkably  similar  for  SLD 
regressions  using  any  of  the  MLD  methodologies.  In  particular,  slope 
values  are  almost  constant  and  negative  in  the  Marmara,  Black  and 
Azov  Seas.  For  the  western  Black  Sea  the  slope  is  nearly  zero  while  the 
intercept  is  positive,  indicating  that  SLD  is  deeper  than  MLD  and  that 
the  difference  varies  little  in  magnitude.  The  Aegean  Sea  has  a 
different  slope  for  MLD(L)  that  may  be  related  to  its  closer  relation¬ 
ship  with  SLD  during  the  spring  in  that  area  (Fig.  5b). 

There  appears  to  be  a  salinity  contribution  to  MLD  in  the  Aegean 
Sea  and  the  northwestern  portion  of  the  Mediterranean  Sea  that  is  in 
our  study  region.  For  the  intercept  associated  with  the  density  based 
MLD  definitions  (MLD(ot)  and  MLD(crc)),  the  Aegean  Sea  tends  to  be 
positive  while  it  is  negative  in  the  Mediterranean  Sea  (Fig.  9a,b). 


Jul  (1784  points) 
2q  —  CoiTelation=  -0.06 


25 


MLD(ap  (m) 


Fig.  7.  Scatter  diagram  of  SLD  versus  MLD(aT)  from  the  January  gridded  field.  The  solid  Fig.  8.  Scatter  diagram  of  SLD  versus  MLD(aT)  from  the  July  gridded  field.  The  solid  line 

line  represents  the  linear  fit  with  a  correlation  of  0.74.  represents  the  linear  fit  with  a  correlation  of  —  0.06. 
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Fig.  9.  Spatial  variation  of  correlation  coefficient  between  SLD  and  for  each  one  of  MLDs  based  on  various  definitions  calculated  over  the  seasonal  cycle.  Similarly,  intercept  and  slope 
values  of  the  least  squares  line  are  computed  over  the  seasonal  cycles  as  further  described  in  the  text. 


Positive  (negative)  indicates  that  SLD  is  deeper  (shallower)  than  MLD. 
For  the  T-only  MLD  definitions  (MLD(r)  and  MLD(L)),  the  northern 
part  of  the  Aegean  Sea  has  a  more  positive  intercept  and  the 
Mediterranean  Sea  is  now  positive  (Fig.  9c, d). 

7.  Summary  and  conclusions 

This  paper  examines  the  evolution  over  the  annual  cycle  of  the  SLD 
relative  to  the  MLD  for  the  Aegean,  Marmara,  Black,  and  Azov  Seas. 
With  historical  quality-controlled  ocean  T  and  S  profiles,  a  monthly 
SLD  climatology  is  constructed  on  a  regular  025°  x  0.25°  grid  using  the 
kriging  methodology.  A  four  step  pre-processing  procedure  designed 
to  reduce  noise  and  the  effects  of  sampling  irregularities  is  applied 
prior  to  forming  the  climatology.  Monthly  SLD  fields  are  then 
compared  relative  to  the  much  more  widely  studied  MLD  as  computed 
using  four  methods  from  the  recent  scientific  literature. 

Very  deep  SLD  occurs  in  the  Aegean  and  Black  Seas  for  the  months 
of  January,  February,  and  March.  This  is  due  to  the  tendency  of  these 
regions  to  have  cooler,  fresher  water  above  warmer,  saltier  water. 
When  the  near  surface  water  becomes  cooler  than  the  deep  water 


below  the  permanent  pycnocline,  the  SLD  drops  to  the  bottom  where 
sound  speed  is  largest.  When  this  occurs,  the  MLD  remains  at  the 
seasonal  pycnocline.  This  is  independent  of  the  methods  used  to 
estimate  MLD,  though  the  MLD(L)  is  generally  shallower  than  MLD 
(cjt)  and  MLD(o'c)  and  has  a  larger  bias.  For  the  rest  of  the  year,  SLD  is 
considerably  shallower  and  closer  to  the  MLD.  During  the  spring  when 
new  stratification  occurs  in  the  Aegean  Sea,  SLD  and  MLD(L)  represent 
the  recent  and  shallow  stratification  while  the  MLD(ot),  MLD(oc),  and 
MLD  (I)  remain  deeper.  Long  range  upper  ocean  acoustic  commu¬ 
nication  would  have  more  favorable  conditions  in  the  months  of 
January,  February,  and  March  in  the  Aegean  and  Black  Seas. 

A  linear  statistical  relationship  is  established  based  on  statistically 
significant  correlation  values  between  SLD  and  MLD.  These  relation¬ 
ships  are  useful  for  (at  least)  two  practical  purposes:  (1 )  to  determine 
the  SLD  from  T-only  profiles  with  the  use  of  MLD(r)  and/or  MLD(L) 
and  the  linear  relationship  and  (2)  when  an  observation  profile  is 
deep  enough  to  estimate  the  MLD  but  potentially  not  deep  enough  to 
estimate  SLD.  Based  on  computation  of  the  error  (c  in  Eq.  (2))  a  perfect 
linear  relationship  (i.e.  correlation  is  ±  1)  results  in  very  small  errors 
(typically <5  m).  In  regions  where  correlations  are  relatively  lower 
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(e.g.  0.4-0.6),  errors  become  larger  as  expected.  In  the  regions  of  the 
Sea  of  Azov  and  mid-western  Black  Sea,  SLD  and  MLD  have  smaller 
magnitudes  (Fig.  5a)  resulting  in  smaller  errors,  within  8  m,  for  all 
seasons  except  summer. 

During  the  winter  when  SLD  is  very  deep,  MLD(L)  has  the  largest 
bias  because  it  tends  to  be  a  shallower  MLD  estimate.  During  the  rest 
of  the  year  MLD(L)  is  closer  to  SLD,  but  the  correlation  coefficient  is 
still  smaller  than  that  for  the  other  MLD  methods.  For  this  reason,  MLD 
(1)  may  be  the  best  predictor  of  SLD,  when  S  is  unavailable,  using  the 
regression  analysis  presented  in  this  paper. 
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Appendix  A 

Four  steps  are  applied  to  the  SLD  and  MLD  data  before  kriging 
interpolation  is  performed.  These  steps  are  designed  to  reduce  noise 
and  even  out  the  relative  weight  of  the  observations. 

A.l  Step  I  Z-score  clipping 

To  minimize  errors,  we  remove  SLD  and  MLD  outliers  statistically. 
A  common  approach  is  to  remove  data  values  that  exist  outside  of  a 
standard  deviation  range.  The  problem  with  this  method  is  that  the 
standard  deviation  computation  itself  is  contaminated  by  the  outliers. 
To  correct  this  problem  we  use  in  place  of  the  standard  deviation  and 
mean  the  corresponding  biweight  estimate  for  the  variance  (scale) 
and  mean  (location),  respectively  (Lanzante,  1996).  The  biweight 
methodology  is  a  robust  statistical  measure  that  accounts  for  non¬ 
normal  distributions  of  data  as  a  result  of  gross  errors  or  outliers.  All 
computation  of  the  mean  and  standard  deviation  are  therefore 
replaced  with  the  biweight  mean  and  standard  deviation  estimates. 
The  removal  of  outlier  data  is  then  achieved  using  the  Z-score  statistic 
described  by  Lanzante  (1996)  and  adopted  by  others  (Zou  and  Zeng, 
2006;  Carrier  et  al.,  2007). 

The  Z-score  is  a  normalization  of  the  data  using  the  biweight  mean 
and  standard  deviation.  For  a  sample  of  n  observations  (x/,  i  =  1,2,3,  -  n) 
the  Z-score  is  given  by 

7.  ^ 


where  Xbi  is  the  biweight  mean  and  Obi  is  the  biweight  standard 
deviation.  Data  are  removed  where  the  Z-score  is  greater  than  4, 
indicating  a  deviation  of  4  biweight  standard  deviations  from  the 
biweight  mean. 

A.2.  Step  2.  Super-observation 

Climatological  super-observations  represent  the  median  value  and 
location  of  a  data  cluster  from  all  years  for  a  given  month  that  are 
within  a  1/12°  circle.  A  set  of  super-observations  are  computed  by 
iterating  through  all  the  original  observations  and  replacing  them 
with  super-observations.  This  is  done  as  follows: 

The  procedure  is  an  iterative  process  for  each  observation 
(i  =  1,2,3,. ..n).  For  the  ith  observation,  all  other  observations  within 
a  1/12°  circle  are  collected,  and  the  median  value  and  median  location 
of  these  observations  (including  the  ith  observation)  are  computed. 
The  ith  observation  and  all  other  observation  within  the  1/12°  circle 


are  then  removed  and  replaced  with  the  computed  median  value  and 
location.  The  discarded  observations  are  not  available  for  the  rest  of 
the  process.  If  there  are  no  other  observations  within  a  1/12°  circle, 
the  ith  observation  is  kept  as  is.  After  one  pass  through  the  obser¬ 
vations,  a  new  set  of  super-observations  exists.  Because  the  median 
locations  of  the  1/12°  circles  shift  relative  to  the  location  of  the  ith 
observation  location,  some  of  the  new  super-observations  are  closer 
together  than  1/12°.  For  this  reason,  additional  iterative  passes  are 
conducted  until  no  further  data  reductions  are  achieved.  This  occurs  in 
3  or  4  passes  with  50%  to  70%  reduction  in  observations,  with  90%  of 
the  data  reduction  occurring  in  the  first  pass. 

A.3.  Step  3.  Median  filter 

The  median  filter  (Tukey,  1977)  is  widely  used  in  image  processing 
for  its  salt  and  pepper  noise  reduction  and  edge  preserving  qualities 
(Maragos  and  Schafer,  1987).  A  constant  time  application  (Perreault 
and  Hebert,  2007)  is  adopted  here  but  for  irregular  data.  Our  median 
filter  iterates  through  each  (super)  observation,  replacing  its  value 
with  the  median  of  the  surrounding  nearest  seven  observations.  The 
location  of  the  observation  is  not  altered.  Seven  observations  were 
chosen  because  for  a  normal  distribution,  sampling  theory  margin  of 
error  is 

F-y  ^ 

where  Zq;/2  is  the  value  for  the  a/2  percent  area  in  the  right  tail  of  the 
normal  distribution,  a  is  the  standard  deviation,  and  n  is  the  number 
of  observations.  This  is  a  two-sided  hypothesis  test  such  that  a  is  the 
percent  probability  that  the  errors  are  less  than  E.  Only  7  data  points 
are  needed  to  achieve  a  95%  confidence  for  a  margin  of  error  that  is 
equal  to  0.75a: 


Since  multiple  applications  of  the  median  filter  are  recommended 
(Gallagher  and  Wise,  1981),  the  median  filter  is  applied  three  times. 
Further  changes  are  generally  small  for  additional  applications. 

A.4.  Step  4:  Smoothing 

While  the  median  filter  removes  salt  and  pepper  noise  and 
preserves  edges,  it  does  not  smooth.  In  order  to  smooth  the  data,  a 
mean  filter  is  applied  to  1/2°  regions  around  each  observation. 
Applied  in  the  same  way  as  the  median  filter,  three  applications 
reduce  noise  while  smoothing  the  horizontal  gradients. 
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